Stripes, Zigzags, and Slow Dynamics in Buckled Hard Spheres 
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We study the analogy between buckled colloidal monolayers and the triangular-lattice Ising an- 
tiferromagnet. We calculate free volume-induced Ising interactions, show how lattice deformations 
favor zigzag stripes that partially remove the Ising model ground-state degeneracy, and identify the 
Martensitic mechanism prohibiting perfect stripes. Slowly inflating the spheres yields jamming as 
well as logarithmically slow relaxation reminiscent of the glassy dynamics observed experimentally. 
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Geometric frustration, manifested for example in the 
anti-ferromagnetic (AF) Ising model on a triangular lat- 
tice [l| , occurs whenever local interaction energies cannot 
be simultaneously minimized. It gives rise to highly de- 
generate ground states, unusual phases of matter [2j, and 
possibly slow or glassy dynamics [3j, whose properties 
even after decades of research are not fully understood. 
Here we present a theoretical study of a colloidal system, 
one of a class of artificial frustrated systems in which the 
state of each constituent can be directly visualized [J], 
that provides new insight into the microstructure of frus- 
trated systems and its connection with their dynamics. 

We study hard spheres confined between parallel 
plates. For plate separation slightly greater than the 
sphere diameter and at sufficiently high sphere density, 
the spheres buckle upward or downward [a, ly, |2] from 
their lower-density positions on a hexagonal lattice. This 
buckling gives rise to a choice of two states for each 
sphere, analogous to the two states of Ising-model spins 
[a, [8|. The tendency of spheres to maximize free vol- 
ume introduces an effective repulsive interaction that fa- 
vors configurations with neighboring spheres in opposite 
states just as the AF Ising interaction favors opposite 
states of neighboring spins. As in the triangular-lattice 
AF Ising model, frustration arises because it is impossible 
to arrange the three particles on any triangular plaquette 
such that all pairs of neighboring particles are in opposite 
states. In the AF Ising model on a rigid lattice, there is 
an extensive number of ground-state configurations (im- 
plying an extensive ground-state entropy) in which neigh- 
boring spins on two of the three bonds on each plaquette 
are in opposite states. Given the analogy just drawn be- 
tween our colloidal system and the AF Ising model, it is 
reasonable to conjecture that the colloidal system might 
exhibit ground-state degeneracies and dynamics similar 
to those of the rigid-lattice AF Ising model. Recent ex- 
periments in diameter-tunable-microgel systems revealed 
sub-extensive ground-state entropy and glassy dynamics 
|9j|. Our theoretical study will address the differences 
between the colloidal system and the rigid lattice Ising 
model and make some conjectures about a likely closer 
analogy between the colloidal system and the AF Ising 
model on a compressible lattice [lOj . 



We show that the short-ranged AF behavior in this 
hard-sphere system may be explained by a simple geo- 
metrical model relating it to the nearest-neighbor Ising 
model. However, the out-of-plane buckling induces local 
in-plane lattice distortions that, as in the elastic Ising 
model, partially remove the Ising ground-state degen- 
eracy and select configurations with zigzagging stripes 
of up and down spheres. This 'ground-state' lacks the 
local zero-energy modes found in the Ising model on a 
rigid triangular lattice, and as a result, the colloidal sys- 
tem exhibits dynamics that are qualitatively slower than 
those of the Ising model. Moreover, stripes require global 
deformations that are incompatible with the system's 
boundary conditions; consequently they break up into a 
Martensite ll| and form a new partially disordered and 
highly degenerate 'ground-state'. 

The free energy of our hard-sphere system is dictated 
by its phase-space volume, which is a collective function 
of all the particles in the system; hence this system is not 
exactly equivalent to an Ising model with pairwise addi- 
tive interactions. Nonetheless, we may compare our sys- 
tem to the nearest-neighbor Ising model on the triangular 
lattice and ask what is the strength of AF interactions 
that best describes the hard-sphere system. 

A 'microscopic' state is specified by the positions 
{xi,yi,Zi} of all particles 1 < i < N. We coarse-grain 
these states into Ising-like configurations specified by {si} 
with Si = sign(zi) (x,y are the coordinates in the plane 
of the confining walls, z is perpendicular to the walls, 
and z — is at the middle of the cell) . The probability 
of a particular configuration {si} is equal to the 3iV- 
dimensional integral V ({si}) over all states belonging to 
that configuration, divided by the total phase-space vol- 
ume Vtot °f a U configurations: 



Pns({si}) = V({si})/V tot . 



(1) 



We would like to equate this to the probability of finding 
the corresponding configuration in the Ising model, 



P l ({ Sl }) = exp 



(pj 






/z, 



(2) 



where /3 = 1/fcgT is the inverse temperature, J is the in- 
teraction strength, the sum runs over all nearest neighbor- 
pairs, and Z is the canonical partition function. 




FIG. 1: (Color online) Free volume model: A) Top view. Con- 
tributions to free volume originate from motion along axes to 
each of the neighbors. B) Side view. Up particle surrounded 
by down particles has more free volume (lower energy) than a 
down particle. C) Surrounding neighbors touch wall and are 
separated by 2L, central particle is confined to the vertical 
plane. Large circles are volumes excluded by neighbors, hor- 
izontal bands are volumes excluded by the walls, remaining 
white region is the central particle's free volume, divided at 
the middle of the cell height (dotted line) into V+ and V- . 



Unlike the commonly used cell-model, which approxi- 
mates the phase-space volume of a system as a product 
of single-particle free volumes, our model equates Phs 
and Pi by assuming that V ({si}) is a product of con- 
tributions from all nearest-neighbor 'bonds': 1^({sj}) = 
Y[ v (siSj] h, d, L) , where the pair contribution v depends 
on SiSj and on the wall separation h, the sphere diameter 
d, and the in-plane number density, which we character- 
ize by the spacing L of the underlying triangular lattice. 
We evaluate v(siSj) in a quasi-one-dimensional approxi- 
mation by allowing particles i and j to move only in the 
vertical plane passing through the axis connecting their 
lattice positions (see Fig. Q3A.). We consider a particle, 
which we call the central particle, and its two neighbors 
along one of the principal lattice directions (see Fig.[lj3). 
If the two neighbors are in opposite Ising states, the free 
volumes resulting from the central one being up or down 
are equal by symmetry. When the two neighbors are in 
the same state (down, without loss of generality), the 
central particle has more free volume (V+) when it is up 
than when it is down (V-). We calculate V+ and V- from 
the geometrical setting of Fig. [Tp and equate the ratio 
of the probabilities of finding the two configurations in 
Fig. QJ3 for hard spheres to that of the Ising model: 



V+ _ exp(-/3P + ) 
V- exp(-/3£L) 



exp(-4/?J). 



(3) 



From this, we deduce that the hard-sphere system corre- 
sponds to an Ising model with an effective AF interaction, 
f3J cS (d,h,L) =-hx(y+/V-)/A < 0. 

For given geometrical parameters d, h, L, we evalu- 
ate V+ and V- and determine the effective interaction 
strength /3J e ff. We then use the exact solution of the 
Ising model [jj| to calculate the average number (Nf) of 
frustrated neighbors per particle (we refer to SiSj = — 1 
as satisfied and to s^Sj = 1 as frustrated), which pro- 
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FIG. 2: (Color online) Average number of frustrated bonds 
per particle vs sphere diameter and cell height. Free-volume 
model (lines) agrees with Monte-Carlo simulations (symbols) . 



vides a measure of short-range AF order. (Nf) = 2 is 
the value in the ground state, and (Nf) = 3 corresponds 
to a random configuration. Figure[2]shows the agreement 
between our simple model and three-dimensional Monte- 
Carlo (MC) simulations. Simulations included N = 1600 
spheres with steps consisting of small displacements of 
single spheres as well as area-preserving box deforma- 
tions, in which the angle or aspect ratio of the simu- 
lated parallelogram was allowed to change [7| . To probe 
cases with d > L, we start with striped configurations 
that can accommodate the maximal sphere diameter by 
lattice deformation (see below), wait 4xl0 5 steps, and 
then average (Nf) over additional 4xl0 5 steps. We plot 
results only of cases for which d was large enough for 
the system to have long-range six-fold orientational or- 
der * 6 = (<exp(i69jk)) > 0.5 (9jk is the angle the bond 
between j and k forms with an arbitrary axis, and the 
average is over all nearest-neighbors pairs [12|]). Small 
spheres (d < L) in a wide cell (h/L = 1.5) are weakly 
confined, and the approximation that the surrounding 
spheres in Fig. [T]C touch the walls fails, giving rise to 
small differences between the model and simulations. 

For large spheres (d> L), simulations remain jammed 
in striped configurations and do not increase the value 
of (Nf) beyond the initial value of 2, even though they 
are expected to do so from the free volume considera- 
tions incorporated in the model. To further explore this 
jamming, we conducted MC simulations that started at 
a disordered configuration with da — L and then ordered 
as the sphere diameter was gradually increased to some 
larger value d. The spheres were initially on a triangu- 
lar lattice in the xy plane with each sphere randomly 
touching either the top or bottom wall. To speed the 
simulations, we considered random jumps in the z direc- 
tion between touching either walls, while keeping the xy 
displacements continuous. During the swelling process, 
once every MC step the diameter of all spheres was in- 
creased to the maximal value allowable without overlaps. 
Figure [3] shows results of simulations with wall separation 
h/L = 1.3. For d/L = 1.005, the free- volume model pre- 
dicts (Nf) = 2.12, and the simulation indeed slowly equi- 
librates to that value by ~10 5 MC steps. For d/L = 1.01, 
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FIG. 3: (Color online) A) Sphere diameter, and, B) average 
number of frustrated bonds per particles following swelling to 
different sphere diameters. Normalized cell height h/L — 1.3. 
Inset: same for Ising model following quench to T = 0. 




FIG. 4: (Color online) Final configurations following swelling 
to d/L = 1.01 at h/L = 1.3. A) Deformable box. B-C) Rigid 
box. System has N = 1600(A-B), 6400(C) spheres. Spheres 
touching top/bottom wall are dark/bright. Simulation boxes 
are deformed parallelograms with periodic boundary condi- 
tions. For ease of presentation we copy the simulated region 
and plot a rectangular region of the periodic system. 



the system's relaxation to the value of (Nf) = 2 pre- 
dicted by the model includes a logarithmically slow de- 
cay to (Nf) rs 2.1 over a time scale of 10 7 MC steps, 
followed by a sharp jump to the equilibrium state. For 
d/L = 1.015, although the system is expected to be at a 
state with (Nf) = 2, it gets jammed during the swelling 
process at a state with (Nf) = 2.35 and does not leave 
it over the time scales investigated here. Note that in 
Fig. [3] we plot a single realization for each case, however 
we observed similar behavior when repeating the simu- 
lations with multiple realizations. Neither jamming nor 
logarithmically slow relaxation occur in the Ising model 
even when quenched to zero temperature (see inset). 

Densely-packed spheres exhibit slower dynamics than 
low-temperature Ising spins on a rigid lattice because the 
morphology of the maximally-packed hard-sphere config- 
urations differ from those of the Ising ground-state. Un- 
like the highly disordered Ising ground-state [1| , the hard- 
sphere 'ground-state' consists of parallel zigzag stripes 
(Fig. HJA.). In the Ising model, each triangular plaqutte 
has one bond frustrated and two satisfied. Although one 
third of the bonds in the system are frustrated, and the 
average number of frustrated neighbors per particle is 
(Nf) = 2, not all particles have exactly two frustrated 
neighbors. By considering the six triangles surrounding a 
certain spin in the lattice, Fig.O'V shows the five possible 
ways (up to rotations and spin inversions) to align them 
such that each triangle will have a single frustrated bond. 
The central spin may have Nf = 0, 1, 2, or 3 frustrated 
neighbors. This leads not only to disorder but also to 
fast relaxation dynamics since spins with Nf — 3 are free 
to flip without an energetic cost. For close-packed buck- 
led spheres, each triplet of spheres in contact defines an 
equilateral triangle with sides d. As in the Ising model, 
one of the three spheres is up (or down) and two down 
(or up), thus tilting this equilateral triangle with respect 
to the horizontal plane. When projected onto the plane, 





FIG. 5: (Color online) Tiling with: A) equilateral triangles 
for the Ising ground-state, and, B) isosceles triangles for close- 
packed buckled hard spheres. The large angle (3 is blackened. 



the tilted equilateral triangle is deformed to an isosceles 
triangle with one long side d along the frustrated bond 
and two shorter sides x = yfd 2 — (h — d) 2 < d along 
the satisfied bonds. Each of these isosceles triangles has 
two small angles a = cos -1 (^-) < -| and a large an- 
gle j3 = 7r — 2a > ?. Now, close-packed configurations 
for the buckled spheres are equivalent to tiling the plane 
with these isosceles triangles. To completely cover the 
plane, the angles of the six triangles meeting at each ver- 
tex must sum to 27r. Figure [SB demonstrates that for 
N f = 0, 1, 3 the angles sum to 6/3 > 27r, 2a+4/3 > 2tt, and 
6a < 27T, respectively, and thus that the triangles can- 
not fit together: for the two configurations with Nf = 2 
the angles sum to 4a + 2/3 = 2-7T, enabling a perfect tiling 
corresponding to the maximal-density close-packed state. 
The slow dynamics observed for large sphere diameters 
result from the lower degeneracy of these zigzagged stripe 
configurations compared to the Ising ground-state. More 
importantly, the close-packed states with Nf = 2 do not 
have the free particles with Nf — 3 that are crucial for 
the low-temperature dynamics in the Ising model [1 31 ] - 
Here, many spheres need to cooperatively rearrange in 
order for the system to find configurations that maximize 
the free volume. Spheres swollen to a very large diame- 
ter (d/L = 1.015 here) hardly move vertically to change 
their Ising-configuration because the neighboring spheres 
do not have enough room to rearrange in the horizontal 



directions and to accommodate the lattice deformations 
required to achieve optimal packing. 

When the spheres swell slowly enough, they find a 
configuration that maximizes free- volume by each sphere 
having exactly two frustrated neighbors. Such configu- 
rations consist of parallel zigzagging stripes (Fig. HJ^). 
Stripes run only along two of the three principal lat- 
tice directions, hence the local distortions are non- 
isotropic and require a macroscopic deformation of the 
system. This is possible in the simulations described 
above in which the shape of the simulation's bounding- 
box changes dynamically [7|. However, experimentally 
the spheres form crystalline domains separated by grain 
boundaries [12j |. which may be better described theoret- 
ically by rigid boundary conditions. Then, the local ten- 
dency for zigzag stripes is incompatible with the rigid 
boundary conditions. The tiling rules for the isosce- 
les triangles induce local deformations along two of the 
three principal lattice directions, and for the system 
to be globally isotropic it must break up into domains 
with stripes running along different directions. We sus- 
pect that this is the mechanism leading to the bro- 
ken stripes seen experimentally _[9J and we indeed ob- 
served such Martensitic states [ll] when repeating the 
swelling simulations without allowing the simulation box 
to deform 14J. For instance, the case of h/L = 1.3, 
d/L = 1.01 relaxes in the deformable box simulations 
to the zigzagged striped state with Nf = 2, whereas in 
a rigid box, the average number of frustrated neighbors 
relaxes to (Nf) — 2.05 (Fig. [3)3), and the final configura- 
tion (Fig. |4}3) consists of broken stripes. We saw similar 
structures (Fig. HP) and relaxation to the same value of 
(N f ) for N = 100,400,1600,6400. It would be interest- 
ing to test whether the size of these domains scales as 
the square root of the system size, as was found in other 
Martensites 11|. and whether subsequently (Nf) slowly 
goes to 2 in the large- N limit. 

The maximal sphere diameter possible in a zigzag con- 
figuration is equal to that of straight stripes, and the 
free-volume-cell approximation does not distinguish be- 
tween the two Nf — 2 configurations corresponding to a 
straight segment of a stripe and to a bend in the stripes. 
However simulations and experiments seem to indicate a 
possible preference for straight stripes over zigzags. It is 
unclear if the observed zigzag patterns represent equiva- 
lence between straight and zigzagged stripes, or whether 
the system falls into zigzagged configurations due to ki- 
netic reasons. It would be interesting to go beyond 
the mean-field description, as was done when comparing 
the face-centred cubic and the random hexagonal-close- 
packed structure of hard spheres in three dimensions 15]. 

The relief of frustration by lattice deformation resem- 
bles the elastic Ising model [10J] , which when analyzed ex- 
actly at the microscopic level yields by our isosceles tiling 
scheme a zigzag-stripe ground-state. It would be interest- 
ing to further investigate the finite temperature behavior 



of that model, as well as other models with zigzag-stripe 
ground-states [16J |. 
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